class: center, middle, inverse, title-slide .title[ # Advanced quantitative methods øvelsestime uge 8 ] .author[ ### Rasmus Klokker ] --- class: center, middle <style type="text/css"> .strike { text-decoration: line-through } </style> <style type="text/css"> .watch-out { background-color: lightpink; border: 3px solid red; font-weight: bold; } .scroll-output { height: 90%; overflow-y: scroll; } .scroll-output-2 { overflow: scroll !important; white-space: nowrap; width: 100%; height: 100%; } .scroll-box-14 { height:14em; overflow-y: scroll; } .scroll-box-14-2 { height:14em; width:14em; overflow-y: scroll; overflow-x: scroll; } .pull-left-1 { float: left; width: 36%; } .pull-right-2 { float: right; width: 60%; } .right-column{ padding-top: 0; } .center2 { margin: 0; position: absolute; top: 50%; -ms-transform: translate(-50%, -50%); transform: translate(-50%, -50%); } </style> #Spørgsmål til undervisningen? --- #multipel regression Nu kan vi kontrollere for confounders!! --- #omitted variable bias Hvad sker der når vi ikke kontrollerer for confounders? *Omitted variable bias*..Det er hvad der sker!! <img src="data:image/png;base64,#6z3g0q.jpg" width="600px" height="400px" /> --- Når der er variable som er blevet "udeladt" fra vores regression, som skulle have været kontrolleret for, så får vi "omitted variable bias"(OVB) Omitted variable bias minder om selection bias i den forstand, at vi "får for meget med" når vi udelader variable, vi ikke skulle have udeladt. I forhold til OVB består bias i at vi også får `\(sammenhæng~mellem~confounder~og~treatment \times samenhæng~mellem~confounder~og~outcome\)` med i købet når vi får en regressionskoefficient. Hvis vi har OVB bliver vores regressionskoeffiecient derfor: `\(\scriptsize \beta_{treatment}=Kausal~effekt_{treatment}+\underbrace{(sammenhæng~mellem~confounder~og~treatment \times samenhæng~mellem~confounder~og~outcome)}_\textrm{omitted variable bias}\)` Det vi var interesserede i at få var jo `\(Kausal~effekt_{treatment}\)`, men indtil vi har kontrolleret for alle confounders så kommer det ikke til at ske. --- #OVB eksempel Lad os tage en tur tilbage til ESS og til variablen "psppipla"(political efficacy). Her vil vi gerne se på effekten af udannelsesniveau("eduyrs") på political efficacy ```r ESS <- read_spss("C:/Users/B059064/Desktop/PHD/teaching/advanced stats/ESS9e03_1.sav") %>% filter(cntry == "DK") %>% # Keep only the Danish data, # Keep only a minimum set of variables we need today, select(idno, pspwght, gndr, eduyrs, agea, psppsgva, trstlgl, trstplc, pdwrk,health, psppipla) %>% drop_na() # Delete cases with missing values. ESS <- ESS %>% mutate( psppsgva = zap_labels(psppsgva), # Make numeric eduyrs = case_when( # Censor years of education at 9 & 21 years. eduyrs > 21 ~ 21, eduyrs < 9 ~ 9, TRUE ~ as.numeric(eduyrs)), gndr = as_factor(gndr)) # Make factor ``` --- ```r mod1 <- lm(scale(psppipla) ~ scale(eduyrs), data = ESS) htmlreg(list(mod1)) ``` <table class="texreg" style="margin: 10px auto;border-collapse: collapse;border-spacing: 0px;caption-side: bottom;color: #000000;border-top: 2px solid #000000;"> <caption>Statistical models</caption> <thead> <tr> <th style="padding-left: 5px;padding-right: 5px;"> </th> <th style="padding-left: 5px;padding-right: 5px;">Model 1</th> </tr> </thead> <tbody> <tr style="border-top: 1px solid #000000;"> <td style="padding-left: 5px;padding-right: 5px;">(Intercept)</td> <td style="padding-left: 5px;padding-right: 5px;">-0.00</td> </tr> <tr> <td style="padding-left: 5px;padding-right: 5px;"> </td> <td style="padding-left: 5px;padding-right: 5px;">(0.03)</td> </tr> <tr> <td style="padding-left: 5px;padding-right: 5px;">scale(eduyrs)</td> <td style="padding-left: 5px;padding-right: 5px;">0.22<sup>***</sup></td> </tr> <tr> <td style="padding-left: 5px;padding-right: 5px;"> </td> <td style="padding-left: 5px;padding-right: 5px;">(0.03)</td> </tr> <tr style="border-top: 1px solid #000000;"> <td style="padding-left: 5px;padding-right: 5px;">R<sup>2</sup></td> <td style="padding-left: 5px;padding-right: 5px;">0.05</td> </tr> <tr> <td style="padding-left: 5px;padding-right: 5px;">Adj. R<sup>2</sup></td> <td style="padding-left: 5px;padding-right: 5px;">0.05</td> </tr> <tr style="border-bottom: 2px solid #000000;"> <td style="padding-left: 5px;padding-right: 5px;">Num. obs.</td> <td style="padding-left: 5px;padding-right: 5px;">1489</td> </tr> </tbody> <tfoot> <tr> <td style="font-size: 0.8em;" colspan="2"><sup>***</sup>p < 0.001; <sup>**</sup>p < 0.01; <sup>*</sup>p < 0.05</td> </tr> </tfoot> </table> --- Men hvad nu, hvis vi troede at alder var en confounder? `\(0.216=kausal~effekt+(sammenhæng~mellem~alder~og~udd. \times sammenhæng~mellem~alder~og~pol.~eff.)\)` Hvad hvis troede at korrelationen mellem alder og udannelsesniveau var dobbelt så stor som korrelationen mellem uddannelsesniveau og political efficacy? Og hvad hvis korrelationen mellem alder og political efficacy var lige så stor? `\(\begin{aligned} 0.216=kausal~effekt+0.216 \times 2 \times 0.216 \times 2\\=\\0.216=kausal~effekt+0.432 \times 0.432 \end{aligned}\)` Hvor stor, eller lille, ville den kausale effekt så være? --- Her er det rigtig heldigt for mig, at R også kan løse ligninger, så vi kan finde `\(kausal~effekt\)` ```r f <- function(x){ x+0.432*0.432 - 0.216 } sol <- uniroot(f, interval = c(-1e+08, 1e+08)) sol$root ``` ``` ## [1] 0.029376 ``` Så bliver effekten af uddannelsesniveau næsten 10 gange mindre! --- Men hvor stærk skal vores confounder være, før effekten helt forsvinder?
--- #OVB opsummering Fortæller os meget præcist, hvorfor vores regressionskoeffiecient ændrer sig som den gør når vi kontrollerer for en confounder(hvis altså vi ved hvad `\(kausal~effekt\)`, `\(sammenhæng~mellem~counfounder~og~treatment\)` og `\(sammenhæng~mellem~counfounder~og~outcome\)` er) kan også bruges til at sige noget om, hvad der ville ske med vores regressionskoefficient, hvis vi en confounder var så og så stærkt korreleret med vores treatment og vores outcome. --- # Hvad sker der når man kontrollerer? #Frisch-waugh-lovell Når man kontrollerer for ting i en regression, er ideen at vi gerne vil holde alt andet en X konstant. Det vil sige, vi vil gerne se hvad der sker med Y, når vi ændrer på X og _intet andet_. Resten af verden skal stå stille mens kun X bevæger sig! --- For at sørge for at alt andet end x "står stille" kan vi gøre brug af residualerne. Residualerne er, oftest, den variation i Y som _ikke_ kan forklares af vores uafhængige variable. Med andre ord er residualerne et slags mål for alt det som foregår i Y, som vores uafhængige variable ikke har noget at gøre med eller har nogen indflydelse på. På den måde kan vi sige at residualerne, `\(\widetilde{Y}\)` er den "del" af Y som ikke har nogen sammenhæng med vores uafhængige variable. <img src="data:image/png;base64,#Week_8_slides_ver2_files/figure-html/unnamed-chunk-4-1.png" width="700px" height="300px" /> --- Hvis vi nu anvender den logik på vores X og confounders, så kan vi sige at Residualerne for X, `\(\widetilde{X}\)`, er X som er "renset" for confounders. <img src="data:image/png;base64,#Week_8_slides_ver2_files/figure-html/unnamed-chunk-5-1.png" width="700px" height="300px" /> --- Hvis vi så også har residualerne fra `\(\widetilde{Y}\)`, så har vi fuldstændigt afskåret vores confounder fra fra både X og Y <img src="data:image/png;base64,#Week_8_slides_ver2_files/figure-html/unnamed-chunk-6-1.png" width="700px" height="300px" /> --- # Hvordan sker det i praksis? Alt det her er heldigvis nemt at gøre i R :) Her prøver vi at se på effekten af at have et haft et betalt arbejde i det sidste stykke tid på ens tillid til det juridiske system. Vi forestiller os, at tilliden til det juridiske system måske lider et knæk når man skal stille bekendtskab med det kafkaske mareridt vi kalder dagpenge-systemet. --- Samtidig er vi dog bekymrede for, om uddannelsesniveau kan være en confounder. Måske er der både højere sandsynlighed for at have et betalt arbejde og for at have tillid til det juridiske system, hvis man har en højere uddannelse? <img src="data:image/png;base64,#Week_8_slides_ver2_files/figure-html/unnamed-chunk-7-1.png" width="700px" height="300px" /> --- Først indlæser vi vores data ```r ESS <- read_spss("C:/Users/B059064/Desktop/PHD/teaching/advanced stats/ESS9e03_1.sav") %>% filter(cntry == "DK") %>% # Keep only the Danish data, # Keep only a minimum set of variables we need today, select(idno, pspwght, gndr, eduyrs, agea, psppsgva, trstlgl, trstplc, pdwrk,health) %>% drop_na() # Delete cases with missing values. ESS <- ESS %>% mutate( psppsgva = zap_labels(psppsgva), # Make numeric eduyrs = case_when( # Censor years of education at 9 & 21 years. eduyrs > 21 ~ 21, eduyrs < 9 ~ 9, TRUE ~ as.numeric(eduyrs)), gndr = as_factor(gndr)) # Make factor ``` --- Og så kører vi den første regression ```r mod0 <- lm(trstlgl ~ pdwrk, data = ESS) mod1 <- lm(trstlgl ~ eduyrs, data = ESS) ``` --- Nu kan vi på en måde "splitte" _trstlgl_ op i to: 1. Den del af variansen i _trstlgl_ som kan forklares af "eduyrs", hvilket er de forudsagte værdier, `\(\widehat{trslgl}\)` Og 2. Den del af variansen som ikke kan forklares af "eduyrs", nemlig residualerne, `\(\widetilde{trstlgl}\)`. Med andre ord har vi den del af _trstlgl_ som udelukkende har "at gøre" med "eduyrs", `\(\widehat{trslgl}\)` , og den del af _trstlgl_ som intet har at gøre med "eduyrs", `\(\widetilde{trstlgl}\)`. ```r #vi får residualerne trst.res <- mod1 %>% pluck(residuals) #vi får de forudsagte værdier trst.fitted <- mod1 %>% pluck(fitted.values) # vi laver variable der indeholder residualerne og de forudsagte værdier ESS <- ESS %>% mutate(trst_res = trst.res, trst_fitted= trst.fitted) ``` --- <img src="data:image/png;base64,#Week_8_slides_ver2_files/figure-html/res-plots-1.png" width="800px" height="500px" /> --- Vi gør det samme med den uafhængige variabel ```r mod2 <- lm(pdwrk ~ eduyrs, data = ESS) #edu.res <- (ESS %>% pull(eduyrs)) - (mod2 %>% pluck(fitted.values)) pdwrk.res <- mod2 %>% pluck(residuals) pdwrk.fitted <- mod2 %>% pluck(fitted.values) ESS <- ESS %>% mutate(pdwrk_res = pdwrk.res, pdwrk_fitted= pdwrk.fitted) ``` --- <img src="data:image/png;base64,#Week_8_slides_ver2_files/figure-html/res-plots2-1.png" width="800px" height="500px" /> --- Når vi sætter det hele sammen ser det altså sådan her ud <img src="data:image/png;base64,#Week_8_slides_ver2_files/figure-html/res-plots3-1.png" width="800px" height="500px" /> --- Og så bruger vi vores residualer som variable ```r ESS <- ESS %>% mutate(trst_res = trst.res, pdwrk_res= pdwrk.res) mod.res <- lm(trst_res ~ pdwrk_res, data = ESS) mod.ctrl <- lm(trstlgl ~ pdwrk+eduyrs, data = ESS) htmlreg(list(mod0,mod.res,mod.ctrl)) ``` <table class="texreg" style="margin: 10px auto;border-collapse: collapse;border-spacing: 0px;caption-side: bottom;color: #000000;border-top: 2px solid #000000;"> <caption>Statistical models</caption> <thead> <tr> <th style="padding-left: 5px;padding-right: 5px;"> </th> <th style="padding-left: 5px;padding-right: 5px;">Model 1</th> <th style="padding-left: 5px;padding-right: 5px;">Model 2</th> <th style="padding-left: 5px;padding-right: 5px;">Model 3</th> </tr> </thead> <tbody> <tr style="border-top: 1px solid #000000;"> <td style="padding-left: 5px;padding-right: 5px;">(Intercept)</td> <td style="padding-left: 5px;padding-right: 5px;">7.45<sup>***</sup></td> <td style="padding-left: 5px;padding-right: 5px;">0.00</td> <td style="padding-left: 5px;padding-right: 5px;">6.19<sup>***</sup></td> </tr> <tr> <td style="padding-left: 5px;padding-right: 5px;"> </td> <td style="padding-left: 5px;padding-right: 5px;">(0.08)</td> <td style="padding-left: 5px;padding-right: 5px;">(0.05)</td> <td style="padding-left: 5px;padding-right: 5px;">(0.20)</td> </tr> <tr> <td style="padding-left: 5px;padding-right: 5px;">pdwrk</td> <td style="padding-left: 5px;padding-right: 5px;">0.39<sup>***</sup></td> <td style="padding-left: 5px;padding-right: 5px;"> </td> <td style="padding-left: 5px;padding-right: 5px;">0.18</td> </tr> <tr> <td style="padding-left: 5px;padding-right: 5px;"> </td> <td style="padding-left: 5px;padding-right: 5px;">(0.10)</td> <td style="padding-left: 5px;padding-right: 5px;"> </td> <td style="padding-left: 5px;padding-right: 5px;">(0.11)</td> </tr> <tr> <td style="padding-left: 5px;padding-right: 5px;">pdwrk_res</td> <td style="padding-left: 5px;padding-right: 5px;"> </td> <td style="padding-left: 5px;padding-right: 5px;">0.18</td> <td style="padding-left: 5px;padding-right: 5px;"> </td> </tr> <tr> <td style="padding-left: 5px;padding-right: 5px;"> </td> <td style="padding-left: 5px;padding-right: 5px;"> </td> <td style="padding-left: 5px;padding-right: 5px;">(0.11)</td> <td style="padding-left: 5px;padding-right: 5px;"> </td> </tr> <tr> <td style="padding-left: 5px;padding-right: 5px;">eduyrs</td> <td style="padding-left: 5px;padding-right: 5px;"> </td> <td style="padding-left: 5px;padding-right: 5px;"> </td> <td style="padding-left: 5px;padding-right: 5px;">0.10<sup>***</sup></td> </tr> <tr> <td style="padding-left: 5px;padding-right: 5px;"> </td> <td style="padding-left: 5px;padding-right: 5px;"> </td> <td style="padding-left: 5px;padding-right: 5px;"> </td> <td style="padding-left: 5px;padding-right: 5px;">(0.01)</td> </tr> <tr style="border-top: 1px solid #000000;"> <td style="padding-left: 5px;padding-right: 5px;">R<sup>2</sup></td> <td style="padding-left: 5px;padding-right: 5px;">0.01</td> <td style="padding-left: 5px;padding-right: 5px;">0.00</td> <td style="padding-left: 5px;padding-right: 5px;">0.04</td> </tr> <tr> <td style="padding-left: 5px;padding-right: 5px;">Adj. R<sup>2</sup></td> <td style="padding-left: 5px;padding-right: 5px;">0.01</td> <td style="padding-left: 5px;padding-right: 5px;">0.00</td> <td style="padding-left: 5px;padding-right: 5px;">0.04</td> </tr> <tr style="border-bottom: 2px solid #000000;"> <td style="padding-left: 5px;padding-right: 5px;">Num. obs.</td> <td style="padding-left: 5px;padding-right: 5px;">1489</td> <td style="padding-left: 5px;padding-right: 5px;">1489</td> <td style="padding-left: 5px;padding-right: 5px;">1489</td> </tr> </tbody> <tfoot> <tr> <td style="font-size: 0.8em;" colspan="4"><sup>***</sup>p < 0.001; <sup>**</sup>p < 0.01; <sup>*</sup>p < 0.05</td> </tr> </tfoot> </table> --- # Øvelse I den her øvelse skal vi, ligesom før, se på sammenhængen mellem "pdwrk" og "trstlgl", når vi kontrollerer for "eduyrs". Men vi er bekymrede for, at folk med dårligt helbred måske i lavere grad har et betalt arbejde samtidig med at de måske også har stiftet bekendtskab med de mere uheldigere dele af bureakratiet i sundhedssystemet og derfor har mindre tillid til det juridiske system. Hvis det er tilfældet, bliver vi nødt til også at kontrollere for helbred <img src="data:image/png;base64,#Week_8_slides_ver2_files/figure-html/unnamed-chunk-8-1.png" width="700px" height="300px" /> --- #Del 1/4 I første del skal vi kontrollere for "eduyrs" på samme måde som vi lige har gjort. Du må godt copy-paste fra slides :) --- #Del 2/4 I næste del skal vi kontrollere for helbred Her skal du/i 1) Køre en regression, hvor du bruger _residualerne_ for "trstlgl" fra Del 1 af øvelsen, som du lige har beregnet og forhåbentlig lavet til en variabel i datasættet, som den afhængige variabel og variablen "health" som den uafhængige variabel 2) Træk residualerne ud fra den model du kørte i 1) og lav en variabel i dit datasæt, som indeholder residualerne. skip til næste slide for at se løsningen --- ```r mod3 <- lm(trst_res ~ health, data = ESS) #edu.res <- (ESS %>% pull(eduyrs)) - (mod2 %>% pluck(fitted.values)) trst.res <- mod3 %>% pluck(residuals) ESS <- ESS %>% mutate(trst_res = trst.res) ``` --- #Del 3/4 Nu skal vi beregne residualerne for "pdwrk". ## Del 3.1 Ligesom før skal vi først kontrollere for "eduyrs", og ligesom før kan i bare copy-paste fra slides og lave en variabel i ESS med de residualer. ##Del 3.2 Nu skal vi så kontrollere for "health". Ligesom før skal du/i bruge de residualer du lige har beregnet, i Del 3.1, som afhængig variabel: 1) Kør en regression, hvor du bruger _residualerne_ for "pdwrk", som du lige har beregnet og forhåbentlig lavet til en variabel i datasættet, som den afhængige variabel og variablen "health" som den uafhængige variabel 2) Træk residualerne ud fra den model du kørte i 1) og lav en variabel i dit datasæt, som indeholder residualerne. skip til næste slide for at se løsningen --- ```r mod4 <- lm(pdwrk_res ~ health, data = ESS) #edu.res <- (ESS %>% pull(eduyrs)) - (mod2 %>% pluck(fitted.values)) pdwrk.res <- mod4 %>% pluck(residuals) ESS <- ESS %>% mutate(pdwrk_res = pdwrk.res) ``` --- #Del 4/4 NU kan vi så se hvad der sker når vi både har kontrolleret for "health" og "eduyrs"! 1) kør en regression hvor du anvender residualerne for "trstlgl" som afhængig variabel og residualerne for "pdwrk" som uafhængig variabel. Gem modellen som objekt 2) kør en regression hvor du ikke kontrollerer for noget, dvs. "trstlgl" som afhængig variabel og "pdwrk" som uafhængig variabel, og ikke mere. Gem modellen som et objekt 3) Kør en regression som i 2), men hvor du kontrollerer for "eduyrs" og "health", Det vil sige "trstlgl" som afhængig variabel, og "pdwrk", "eduyrs" og "health" som uafhængige variable 4) Sammelign resultaterne fra de 3 regressioner, evt. i "screenreg". Skip til næste slide for at se løsningen --- ```r mod.naive <- lm(trstlgl ~ pdwrk, data = ESS) mod.res <- lm(trst_res ~ pdwrk_res, data = ESS) mod.ctrl <- lm(trstlgl ~ pdwrk+eduyrs+health, data = ESS) screenreg(list(mod.naive,mod.res,mod.ctrl), digits=4) ``` ``` ## ## ===================================================== ## Model 1 Model 2 Model 3 ## ----------------------------------------------------- ## (Intercept) 7.4505 *** -0.0000 7.1637 *** ## (0.0797) (0.0492) (0.2379) ## pdwrk 0.3876 *** 0.0356 ## (0.1035) (0.1063) ## pdwrk_res 0.0387 ## (0.1063) ## eduyrs 0.0888 *** ## (0.0138) ## health -0.4014 *** ## (0.0562) ## ----------------------------------------------------- ## R^2 0.0093 0.0001 0.0725 ## Adj. R^2 0.0087 -0.0006 0.0706 ## Num. obs. 1489 1489 1489 ## ===================================================== ## *** p < 0.001; ** p < 0.01; * p < 0.05 ``` --- #fortolkning Fortolkningen er af vores regressionskoefficient er heldigvis fuldstændig ligesom før. Nu kan vi dog tilføje en enkelt ting: Nu kan vi sige, at vores regressionskoefficient angiver ændringen der sker i vores outcome, for en enhedsændring í vores treatment, _når vi kontrollerer for vores kontrol variable_ I vores eksempel kan således sige at regressionskoefficienten for "pdwrk" er ændringen i "trstlgl" når vi øger "pdwrk" med en enhed, _kontrolleret for "eduyrs" og "health_ --- ## Kan vi fortolke vores kontrolvariable? Når vi har kontrolleret for vores kontrolvariable, tror vi jo ikke der mere confounding på `\(X_1\)`. Men derfor kan vores kontrolvariable jo godt være udsat for confounding <img src="data:image/png;base64,#Week_8_slides_ver2_files/figure-html/unnamed-chunk-9-1.png" width="700px" height="300px" /> --- #gennemgang af øvelse --- # Vi ses næste uge!! ## Office hours 16.2.08, mandag-torsdag 13-15 rhk@vive.dk rhk@soc.ku.dk